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Abstract 

The distribution of inter-occurrence time between seismic events is a 
quantity of great interest in seismic risk assessment. We evaluate this 
distribution for different models of earthquakes occurrence and follow two 
distinct approaches: The non homogeneous Poissonian and the non Pois- 
sonian one. In all cases we obtain either a power law or a power law 
damped by an exponential factor behaviour. This feature of the distri- 
bution makes impossible any prediction of earthquakes occurrence. Nev- 
ertheless it suggests the interpretation of the earthquake occurrence phe- 
nomenon as due to some non-linear dynamics to be further investigated. 

1 Introduction 

Seismology can be defined as the science of earthquakes and studies mainly the 
physics of the earthquakes sources, the elastic wave propagation and the occur- 
rence of earthquakes in space, time and energy. The investigation of earthquakes 
sources and wave propagation are based on the analysis of the seismograms un- 
der the assumption that linear theory of elasticity holds and are treated in a 
huge amount of literature (see e.g. Ref. |I] and references therein). On the 
other hand the study of earthquake occurrence regards the construction and 
the investigation of seismic catalogues, whose content is generally the time of 
occurrence, the location and the energy of earthquakes. The great interest dedi- 
cated by many researchers to the statistics of these quantities is obviously linked 
to the interest in predicting the time, the location and the energy of the next 
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earthquake. These questions are still rather unclear and we shall discuss some 
aspects at the origin of the question. In particular, we suggest that difficulties 
in prediction are intrinsic to the occurrence mechanism. 

The energy release in a seismic event is generally expressed by the mag- 
nitude, which is proportional to the logarithm of the maximum amplitude of 
the recorded seismic signal. The distribution of magnitudes is described by an 
exponential law usually referred as the Gutenberg-Richter distribution, by the 
names of the researchers who firstly observed this feature of seismic catalogues. 
They found that the frequency of occurrence of earthquakes with magnitude 
greater then M behaves as 

logiV = a-bM (1) 

where a indicates the overall seismicity and b is a scaling parameter which 
has typically values close to one. Fluctuations (up to 30%) of the value of b 
around its typical value are widely observed depending on the catalogue, the 
estimating method and the magnitude range 2,. Regional variation of the b 
value are also observed by many authors (see e. g. 0E])- Such features could 
be an indication that scaling properties of earthquakes are not universal. 

However, it has been suggested that this discrepancy in b values could be 
due to systematic errors in magnitude determination [o1|d]. For this reason, it is 
often studied the distribution of the seismic moment M defined as M = fiAAu 
where fi is the shear modulus, A is the area of the seismic fault involved and Au 
is the slip of the fault due to a single seismic event EJ El 03 ED • This is a more 
physical quantity compared to the magnitude and can be obtained by inverting 
the seismic signals with a standard procedure |10|. The cumulative distribution 
of the seismic moment exhibit a power law behavior followed by a sharp cutoff 
after a M™ ax well represented by a Gamma distribution 

-Mq 

/(M ) = CM- (1+7) e^ (2) 

where C is a constant and 7 is a scaling parameter, having a very stable 
value 7 ~ 0.7 EHIIl. 

The fractal nature of the spatial distribution of earthquakes has been shown 
for the CMT catalogue |121 IT3"| and the fractal dimension Df of the hypocenter 
distribution has been determined. However, the fractal dimension could be not 
sufficient to describe all the scaling properties in a given problem and a spectrum 
of dimensions turns out to be necessary in order to fully characterize the scaling 
behavior. If this happens, the distribution is said to be multifractal |2|.It has 
been shown ^JE3 *^ a * ^ e s P a tial distributions of earthquakes in Kanto region, 
in eastern Japan, in California and in Greece have a multifractal structure. The 
multifractal distribution of earthquakes hypocenters has been also confirmed for 
some Italian regions an d it has been suggested ^JH2S] that the temporal 

changes of Df could be a good precursor parameter for earthquake occurrence 
prediction. 

On the other hand, the rate of occurrence of seismic events in time has been 
widely investigated (see e.g. [22 1221 123| ) because the existence of a charac- 
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Figure 1: The experimental distribution of the waiting times for a) the Califor- 
nian catalogue of earthquakes (De Natale et al., 2003), b)the Colfiorito (Italy) 
sequence(De Natale et al., 2003). 

tcristic time scale could make possible the prediction of the next earthquake. 
Unfortunately this is not the case since the distribution of waiting times between 
earthquakes exhibits a two power law behavior |24l 1251 127>| (figure 1). More pre- 
cisely, the data exhibit a first power law regime with an exponent close to 0.8 
for both catalogues, followed by a second power law with an exponent 1.8 for 
Colfiorito and 2.6 for California. The experimental results suggest that inter- 
arrival times are possible at all time scales, making the prediction of earthquakes 
occurrence extremely difficult. The non-Poissonian behavior of earthquake oc- 
currence is clearly due to the phenomenon of clustering, suggested by the power 
law distribution of the inter-arrival times. In fact, for a Poissonian process, the 
distribution would have an exponential behavior. 

The existence of power law behaviors and the multifractal distribution of 
hypocenters, leads to the interpretation of earthquake as a critical phenomenon 
|27U28ll2l?] . proposing a new point of view for the features of earthquake occur- 
rence. 

Finally in recent years it has been proposed that Self Organized Criticality 
(SOC) 30 could explain the observed self - similar properties of earthquakes 
|27l I28| . in particular could be able to reproduce the earthquake size distri- 
bution. Within this approach, the lithosphere structure derives from the self- 
organization of the earth crust in a continental plate. A field theory can be also 
derived from symmetry and conservation laws to explain the earthquakes size 
distribution and viewing the earthquakes as fluctuation of the elastic energy in 
the crust [3"T] . 

In this paper we focus on the waiting time distribution, we derive the ana- 
lytical expression of the inter-arrival time distribution for some existing models 
and we discuss the scaling behavior of the distributions depending on parame- 
ters. More precisely, in Section 2 we shortly present some of the most commonly 
used models and, in section 3, we derive the analytical expression of the prob- 
ability density function using a non homogeneous Poissonian approach for the 
Omori law and the ETAS (Epidemic Type Aftershock Sequences) model. In 
section 4 we derive the p.d.f. using a non Poissonian approach for the Omori 
law and the Poisson Generalized model. In all the cases there exists a choice of 
the parameters fitting some data set. In the final section we present conclusions 
and considerations for earthquakes occurrence prediction. 

2 The earthquake clustering and the Omori law 

It is widely observed that earthquakes tend to occur in bursts. These bursts may 
start suddenly immediately following a large main event, giving rise to the well 
known main - aftershocks sequences, or may build up and die very gradually in 
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time, generating swarms of events. The most important implication of this kind 
of occurrence is that we cannot assume a Poissonian occurrence of earthquakes, 
where a Poissonian process is characterized by a constant rate of occurrence, 
but rather a clustered one. In his pioneering paper, Omori [32] investigated the 
problem of earthquake occurrence within a single cluster of events and proposed 
that the non-Poissonian behavior of seismic catalogues could be well fitted with 
the Omori law, stating that the number of aftershocks n(t) decays in time as 

n{t) = ^ (3) 

where p is generally very close to 1 ranging from 0.7 to 1.7, c is an initial 
time introduced in order to avoid the divergence at t = and k = n(0)c p is an 
experimental constant. 

A widely used approach to earthquakes clustering is provided by "trigger 
model" jSJJ. This assumes a Poissonian occurrence of main events, whereas the 
occurrence of the "triggered" earthquakes is described in terms of a correlation 
function g(t — ti), where U is the time of occurrence of the i-th event. The 
function g(t — U) describes the correlation of each event occurring at time t 
with all the events occurred at previous times. Thus the rate of occurrence will 
be 

* = M+ 9{t-U) (4) 

i:ti<t 

where a is the Poissonian rate of the main events. Among the trigger mod- 
els a widely used one is the Poisson Generalised model this assumes the 
sequence of events as composed by uncorrelated main events which generate 
clusters of aftershocks distributed as the Pareto power law |33J 

«0> - g (5) 

where £(/?) is the Riemann function and j is the number of events in the 
cluster. This approach has been applied for many areas of the world, as Cali- 
fornia js2], Messina Strait area and Campi Flegrei (Italy) [HSj) m order to 
determine the j3 value, found to be between 2.5 and 4. 

A more appropriate choice of g{t — ti) is provided by the ETAS model |37| . 
which considers the existence of many clusters described by the Omori law. 
The model states that the intensity function (the rate of occurrence) of the 
earthquakes is given by: 

foam-Mo) 

a = m+ E u : : ^ ( 6 ) 



i-U< 



t {t-U + c)i 



where u is again the Poissonian rate , a an experimental constant, Mo is the 
smallest magnitude in the catalogue and Mi is the magnitude of the i-th event. 
The meaning of equation (6) is that each earthquake can generate "its own 
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aftershocks" and that the number of these aftershocks depends exponentially 
on the magnitude of the " main" . In other words the clustering degree varies in 
time, leading to a clustering within the clustering. 

A completely different approach is the fractal one A Poissonian pro- 

cess would fill stochastically all the temporal axis and thus would have a fractal 
dimension equal to 1, whereas a clustered process is characterized by a frac- 
tal dimension less then one depending on the clustering degree. Using the box 
counting method it has been found that the New Hebrides seismicity is clus- 
tered with a fractal dimension ranging between 0.126 and 0.255 25:. Moreover 
other authors .38 found that many catalogues in the world have a multifractal 
distribution of inter-arrival times. This result is in good agreement with the 
predictions of the ETAS model. 



3 The non homogeneous Poissonian approach 

The problem of earthquakes inter-arrival time distribution has never been treated 
from the theoretical point of view. In order to explain the temporal clustering 
properties of seismic events occurrence, the most of efforts were dedicated to 
the study of the rate of occurrence. Nevertheless the waiting time distribution 
is very important in the seismic risk assessment because it is very useful in the 
definition of the probability of the occurrence of next earthquake. In this section 
we derive the distribution of the waiting times for a single cluster following the 
Omori law and for the ETAS model. 

The cumulative distribution of waiting times F(At) can be written as |39j 

F(At) = 1 - F (At) (7) 

where Fo(At) is the probability of observing zero events in At. Since for a 
Poissonian process 

Fo(Ai) = e~^ At (8) 
the probability density function (p.d.f.) is found to be 

which is the well known result for a Poissonian process. 

This approach can be generalized also for processes for which /i is not con- 
stant in time and we shall have a non homogeneous Poissonian process. In this 
case the probability of having zero events in At is given by 

F (Ai) = /o At A(t)dt (10) 

where X(t) is the time dependent rate of occurrence thus inserting the (10) 
into the (7) and the result into the (9) we obtain the waiting times p.d.f. Note 
that the Poissonian behavior is more restrictive then the independence of events, 
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Figure 2: The p.d.f. of the waiting time within a cluster of events for different 
values of the model parameter 



since it is obtained under the assumption that the probability of observing more 
than one event in any small time interval, is negligible. 

As a first application of this approach we shall derive the waiting time dis- 
tribution within a cluster of events. In this case the rate of occurrence X(t) is 
given by the Omori law, thus for p ^ 1 

F (At) = e~iT T^W dt (11) 
and the p.d.f. of the At will be 

/(At) = ke^ (At + c)- p e^ M+ ^- p (12) 

which, except for some constant factors, is a Weibull distribution, i.e. a 
power law damped by a stretched exponential decay. Figure 2 shows the p.d.f. 
for different parameter value: the Omori law exponent p controls the decaying 
exponential factor which is dominant at long times for (At + c) 1 ~ p > -^r 2 when 
p < 1 and at short times for (At + c) 1_p < when p > 1. We observe that 
the p.d.f. for p > 1 are not in agreement with experimental observations. 

In the case p = 1 from equation (11) we obtain 



jma.) - (*±£)- (i 3 , 



and 



/(A t ) = fcc*(^±£)- fc - 1 (14) 

which is a power law and does not exhibit any exponential decay as equation 
(12). 

A more complex formula is obtained if we adopt the ETAS model. In this 
case we consider the existence of many clusters of events as described in section 
2. The rate of occurrence is given by the equation (6). If we take the continuum 
limit, that is 



. . k(M) f l - At k(M) 

hit (*-**+ °y ~* Jo ( t - r + c ) p 

where k(M) = ke a (M l - M ) we will get for p ^ 1 



dr (15) 



Fo(Ai) - e~» At+ (i-X-p) {(^+°) 1 -'[a«(i-p)+c]-c»-»»} (16) 



which provides for the p.d.f. 
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Figure 3: The p.d.f. of the waiting times for the ETAS model and p ^ 1. 



Figure 4: The p.d.f. of the waiting times for the ETAS model in the case p = 1 
for different values of k 



/(At) = fo-^(At + C )-*(A* + ^)] 

e -^At+ Tr ^^ T {(At+c) 1 -"[At(l- P )+c]-c 2 -^ 

Equation (17) is well defined i.e. is a positive quantity, only for p > 1 and 
assumes the shape of a Weibull distribution (figure 3). 
Analogously the case p — 1 gives 

F Q (At) = c ck (At + c )-c* e -A^+fcin2±£i) ( 18 ) 
and for the p.d.f. 

/(At) = [/i + kln^^] c ck (At + c )-ck c -At( M+ ki„^M) (1Q) 

We find again a power law damped by an exponential factor. Note that the 
term in the square brackets is negative for \i < fcln ^ 1 because c + At < 4 At, 
however it is possible to obtain positive values of the p.d.f. for k < 0.36 if we 
set \x = 0.5 (figure 4). 

It is noteworthy that the non homogeneous Poissonian approach does not 
provide a good agreement with experimental data since does not predict the 
two power regime shown in figure 1. This feature could be due to the Poisso- 
nian assumption which assumes negligibly small the probability of two events 
occurring in any small time interval. 



4 The non Poissonian approach 

In this Section we derive the analytical expression of the p.d.f. assuming only 
that the probability of cluster occurrence is independent on the probability of 
earthquake occurrence within a cluster. If we call Q n (At) the probability of 
having n events in a cluster and Pjv(At) the probability of having N clusters in 
At, we will have 

F (At) = P (At)Q n (At)+ 

P (At)[l - Qn(At)] + P N (At)[l - Q n (At)} (20) 

The three terms in equation (20) represent respectively the probability of 
having zero clusters of n events, zero clusters of zero events and N clusters of 
zero events. Firstly we determine the p.d.f. of the At within a single cluster. 
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In this case P/v(Ai) = 1 and Po(Ai) = 0. The number of events j in a time 
interval r for p ^ 1 will be given by 

Assuming the power law distribution (5) for j, we have 

1 " 

Qnm = W) Y { r (22) 

Noticing that in the continuum limit — ► J dt and neglecting the quantity 
c 1_p , we have 

«»<^M(T^)X,A« (T + crW1 ~* W = 

[(At + ^-(l- P ) — ] (23) 
where (5=1— /3(1 — p). Finally we obtain the p.d.f. 

'w=w)(Thy f{At+c) - m -' ] « 24) 

which is a power law well defined only for p < 1. This constraint is due to 
the assumption that c 1_p <C (r + c) 1_p , which implies from (21) that, if p > 1, 
the number of events j would became negative. In the case p = 1 we obtain a 
p.d.f. whose behavior is inconsistent with the experimental data and thus will 
not be reported here. 

Next we apply the non Poissonian approach to the "trigger" model which 
assumes a Poissonian occurrence of clusters and a power law decrease of the 
number of events within the clusters (equation (5)). Under these assumptions 
we have 

mA* 



Observing that for a Poissonian process the total number of clusters is TV = 
p,At and that (p.At)\ in the continuum limit becomes T(l + fiAt), we have 

— T T 71 

p (5)-^[ 1 + raW cw -£'"')] m) 

where r = /j,At. Using equation (22) in order to evaluate Q n (At) and 
neglecting again the quantity c^ 1_p ^ we obtain 



8 



Figure 5: The p.d.f. of the waiting times for the Poisson generalized model at 
fixed values of (3 = 3.5 and p — 0.85 and for different values of 

Figure 6: The p.d.f. of the waiting times for the Poisson generalized model at 
fixed values of (3 = 3.5 and [i = 0.01 and for different values of p. 

[K + , m p) + a^)(^±^-lnr)}} (27) 

where a = S = 1 - /3(1 - p), K = aS, b = + cS 1 ^ ~ and 

(j)(x) = (b + c)~ s — (c + x) s . Equation (27) is a very complex function and 
does not allow any simple fit of experimental data. Moreover the number of 
parameters involved in the function is too high for a stable fit. However we 
notice that it is possible to find some plausible relations among some of the 
parameters. For instance, /i and k can be related since they are both rates of 
occurrence: the first one concerns the cluster occurrence, whereas the second 
one states how many earthquakes occur at the beginning of a given cluster. In 
order to evaluate equation (27), we choose k = 20/i. Any other choice for k 
and \i does not influence the shape of equation (28), but only the level of the 
seismicity, that is the total number of events. Obviously the value of fi and k, 
representing the time scale in the system, implies as a consequence the value of 
the constant c of the Omori law and therefore we choose c = 0.3/fc. This means 
that we have three free parameters /i, (3 and p. By varying these parameters, 
we obtain two possible behaviors: either a two power law regime or a two power 
law regime damped by an exponential decay at high At. 

In figure 5 we show the p.d.f. with fixed (3 and p for different values of [i. 
At lower values of fi, therefore for clusters more sparse in time, we observe the 
two power law regimes behavior, whereas for increasing [i we observe the onset 
of an exponential cut-off at long waiting times. Note that the exponents arc in 
the range 0.2 - 0.5 for the first power law and 1.3 - 1.9 for the second one. 

Figure 6 shows the behavior of p.d.f. at fixed (3 and /i for a range of values of 
p. In this case we observe for decreasing values of the Omori exponent p, i.e. for 
clusters lasting a longer time, the onset of an exponential cut-off at long waiting 
times after the two power law regime. In this case the exponents vary between 
0.4 and 0.7 for the first power law and between 1.0 and 1.3 for the second 
one. Any variation of the parameter values does not change substantially the 
behavior in figures 5 and 6. On the contrary, we will sec that the p.d.f. function 
is more sensitive to combined variations of (3 and p. 

Figure 7 shows for p = 0.75 the onset of an exponential cut-off at long 
waiting times for high values of (3 as observed in figures 5 and 6 (slopes are in 
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Figure 7: The p.d.f. of the waiting times for the Poisson generalized model at 
fixed values of p = 0.75 and fj, = 0.01 and for different values of (3. 

Figure 8: The p.d.f. of the waiting times for the Poisson generalized model at 
fixed values of p — 0.95 and /x = 0.01 and for different values of (3. 

the ranges 0.5 - 0.7 and 1.6 - 1.4). On the other hand, for p = 0.95 the two 
power law behavior is substantially insensitive to (3 variations (figure 8). This 
suggests that scaling properties of equation (27) are dominated mainly by p 
than by (3. In this case the power law exponents are 0.4 and 1.0. 

The two power law regime is widely observed for many catalogues in the 
world (figure 1) and generally interpreted as due to catalogue incompleteness. 
Within the Poissonian Generalised approach we find that the two power law be- 
havior is quite robust with respect to parameter changes. Therefore we suggest 
that this feature is an intrinsic property of earthquake occurrence related to the 
P. G. model. 

5 Conclusions 

We evaluate the probability density function of the inter-occurrence time be- 
tween earthquakes following two different approaches. We first assume a non 
homogeneous Poissonian behavior and find for different models of earthquakes 
occurrence always a single power law, eventually followed by an exponential 
decay. 

Next we investigate a non Poissonian approach for different models. The 
obtained p.d.f. has a power law behavior in the case of a single cluster of 
events described by the Omori law. On the contrary, in the case of the Poisson 
Generalised model the p.d.f. exhibit a more complex behavior depending on 
parameters. For all values of p ^ 1 we find consistently a two power law regime. 
This situation, occurring for small fi, corresponds to long waiting times between 
clusters of seismic events, which is the situation more frequently observed in 
nature. Depending on parameters, the value of the exponents are in agreement 
with the experimental data. 

Moreover, for high values of (3, i.e. fast decay in the number of events in 
a single cluster, and a high Poissonian rate n the two power laws are followed 
by an exponential decay. This feature characterizes a weak clustering in the 
distribution of events in time or a frequent cluster occurrence. 

The two power law behavior is observed for many catalogues relatives to 
different areas in the world. This feature, often interpreted as a sign of the 
incompleteness of the catalogue, is here obtained as a specific characteristics of 
the p.d.f. for the Poisson Generalised model. Finally we notice that for all the 
discussed approaches and model the power law behavior implies the absence of a 
characteristic inter-occurrence time and therefore impossibility of any prediction 
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of earthquake occurrence. 

This work is part of the project of the Regional Center of Competence "Anal- 
ysis and Monitoring of the Environmental Risk" supported by The European 
Community on Provision 3.16. 
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